#include <iostream>
#include <stdlib.h>

#include <cuda.h>

#include "libmri.h"

using namespace std;

int main (int argc, char **argv)
{
	int direction = CUFFT_FORWARD;

	if (argc > 1){
		if (!strcmp (argv[1], "forward")) direction = CUFFT_FORWARD;
		if (!strcmp (argv[1], "inverse")) direction = CUFFT_INVERSE;
	}

	cout << "direction = " << direction << " ("
		<< (direction == CUFFT_FORWARD ? "FORWARD" : 
		    (direction == CUFFT_INVERSE ? "INVERSE" : "I DUNNO"))
		<< ")" << endl;

	cuda(SetDevice(0));
	shm_array_t *shm_input = 0;
	char *output = 0;

	if (0 == shm_input){
		if (get_shm_array("input", &shm_input)){
			perror("input");
			return -1;
		}
	}

	if (0 == output){
		output = strdup("output");
	}

	shm_array_t *shm_output = 0;
	if (make_shm_array (output, shm_input, &shm_output)){
		perror(output);
		return -1;
	}

	int xs1 = shm_input->dims[0], xs2 = shm_input->dims[1], batch = 1;
	for (int i = 2; i < shm_input->ndims; ++i){
		batch *= shm_input->dims[i];
	}

	cufftHandle plan;
	plan_cufft2c (plan, xs1, xs2, batch);

	cplxf *x;
	cuda (Malloc (&x, sizeof(cplxf)*xs1*xs2*batch));
	cuda (Memcpy (x, shm_input->data.c, sizeof(cplxf)*xs1*xs2*batch,
	              cudaMemcpyHostToDevice));

	cufft2c_exec (plan, direction, x, xs1, xs2, batch);

	cuda (Memcpy (shm_output->data.c, x, sizeof(cplxf)*xs1*xs2*batch,
	              cudaMemcpyDeviceToHost));
	cuda (Free (x));
	return 0;
}
